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ABSTRACT 

The core of a neutron star contains several species of particles, whose relative equi- 
librium concentrations are determined by the local density. As the star spins down, its 
centrifugal force decreases continuously, and the star contracts. The density of any given 
fluid element increases, changing its chemical equilibrium state. The relaxation towards 
the new equilibrium takes a flnite time, so the matter is not quite in chemical equilibrium, 
and energy is stored that can be released by reactions. For a neutron star core composed 
of neutrons (n), protons (p), and electrons (e), the departure from chemical equilibrium is 
quantifled by the chemical potential difference = /ip + /ie — /Un • A flnite Sfx increases the 
reaction rates and the neutrino emissivity. If large enough {\Sfj,\ > 5A;T), it reduces the net 
cooling rate because some of the stored chemical energy is converted into thermal energy, 
and can even lead to net heating. A simple model (for non-superfluid matter) shows the 
effect of this heating mechanism on the thermal evolution of neutron stars. It is particu- 
larly noticeable for old, rapidly spinning stars with weak magnetic flelds. If the timescale 
for variations of the rotation rate is much longer than the cooling time, a quasi-equilibrium 
state is reached in which heating and cooling balance each other and the temperature is 
completely determined by the current value of P / (or the spin-down power). If only 
modifled Urea reactions are allowed, the predicted quasi-equilibrium x-ray luminosity of 
some millisecond pulsars approaches the upper limits obtained by Banner et al. (1994) 
from ROSAT data. The predicted x-ray luminosity is much lower if direct Urea or other 
fast reactions are allowed. In both cases, the luminosity is probably increased if the stellar 
interior is mostly superfluid. 

Subject headings: dense matter — equation of state — stars: neutron — stars: pulsars: 
general — stars: rotation — x-rays: stars 
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1. INTRODUCTION 

The theoretical study of the thermal evolution of neutron stars and the comparison of 
its predictions to observations of the surface emission of neutron stars of known ages can 
(at least in principle) serve as an interesting and powerful diagnostic of the properties of 
matter at supernuclear densities. The main reason for this is that the dominant cooling 
mechanism of a hot neutron star is the emission of neutrinos and antineutrinos produced 
in the so-called Urea processes (Gamov & Schoenberg 1941): 

n > p + e + i^e, 

p + e > n + Ue- 

The symbols n, p, e, i^e, and i>e denote neutrons, protons, electrons, electron neutrinos, and 
electron antineutrinos. The neutrino emissivity depends sensitively on the microphysics 
of the neutron star core. The reactions can proceed as written ("direct Urea reactions") 
only if the ratio x = rip/n of the proton number density to the total baryon number 
density exceeds a certain threshold value which allows simultaneous energy and momentum 
conservation if all reacting particles are on the surface of their respective Fermi seas (Boguta 
1981; Lattimer et al. 1991). Otherwise, an additional particle participates in the reaction 
and absorbs the additional momentum, without changing its "flavor" ("modifled Urea 
reactions;" Chiu & Salpeter 1964). If the additional particle is another neutron or proton, 
the rate is strongly suppressed, because of large phase space reduction due to degeneracy. 
The suppression is smaller if the momentum can be absorbed by bosons such as pions 
(Bahcall & Wolf 1965) or kaons (Brown et al. 1988). These and other variants of Urea 
reactions are reviewed by Pethick (1992). 

In addition to Urea reactions, other neutrino emitting processes (described in Nomoto 
& Tsuruta 1987) contribute to the cooling, as does the emission of thermal x-ray and 
ultraviolet photons from the stellar surface, which dominates for cooler stars. All cooling 
processes are strongly affected by Cooper pairing (superfluidity or superconductivity) of 
the baryons, which reduces the reaction rates and (at low enough temperatures) the heat 
capacity of the paired particles (e.g.. Wolf 1966; Maxwell 1979; Levenflsh & Yakovlev 
1994a, b). 

Some heating mechanisms may also be present. In isolated neutron stars, the most 
important source of energy that is likely to be tapped is the kinetic energy of rotation, 
which can be converted into internal heat in several ways. If the star does not rotate as 
a solid body, the energy in the differential rotation can be dissipated. This is predicted 
by the standard models of pulsar "glitches," which involve a crustal neutron superfluid 
with pinned vortex lines (Alpar et al. 1984; Shibazaki & Lamb 1989; Umeda et al. 1993; 
Van Riper, Link, & Epstein 1994). Differential rotation and the resulting dissipation 
may also occur in the core of the star, especially if the latter is stably stratifled (Pethick 
1991; Reisenegger & Goldreich 1992) and therefore cannot be spun up by the standard 
Ekman pumping process (Walin 1969; Sakurai 1969). The differential rotation may still 
be impeded by a strong magnetic fleld threading the core, but the existence of a magnetic 
fleld in the core is not well established, and the spin-down process in a strongly stratifled, 
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conducting fluid in the presence of a magnetic fleld has not yet been studied. In addition, 
as the star slows down, its equihbrium shape changes, and the response of the sohd crust 
to this change, through cracking or plastic flow, also generates heat (see Cheng et al. 1992 
and references therein). 

The purpose of this paper is to point out another channel by which a neutron star's 
rotational energy can be converted into heat, and to study its effects. The physical idea 
is as follows. The core of the star contains several kinds of particles, whose relative con- 
centrations adjust by Urea reactions. The equilibrium concentration of particles of each 
kind depends on density. As the star spins down, the centrifugal force decreases continu- 
ously, increasing the star's internal density. This changes the chemical equilibrium state 
throughout the core. Since the relaxation of the relative particle concentrations requires 
a flnite amount of time, the star is never exactly in chemical equilibrium, so energy is 
stored. In addition, the flnite departure from chemical equilibrium modifles the reaction 
rates (Haensel 1992). If the departure from equilibrium is large enough, the net effect 
of the reactions is to increase the thermal energy at the expense of the stored chemical 
energy. 

The dissipation associated with reactions produced by density changes is often called 
"bulk viscosity" or "second viscosity" (Landau & Lifshitz 1959), and has been studied 
as a damping mechanism for oscillations and gravitational-wave instabilities of neutron 
stars (Finzi & Wolf 1968; Sawyer 1989; Cutler, Lindblom, & Splinter 1990; Haensel & 
Schaeffer 1992; Reisenegger & Goldreich 1992; Lai 1994)"^, and as a heating process in 
the collapse of a neutron star to a black hole (Gourgoulhon & Haensel 1993), but not (to 
my knowledge) in the context in which it is discussed in the present paper. I will avoid 
the terms "bulk viscosity" and "second viscosity," since they seem to imply a dissipation 
process that depends locally in time and in space on the divergence of the velocity fleld, 
as suggested, e.g., by equation (15.6) of Landau & Lifshitz (1959). Dissipation by shear 
viscosity does, in fact, depend locally on the shear of the velocity fleld (if the mean free 
path of the particles is much smaller than the macroscopic scales of interest), but "bulk 
viscosity" can be substantially more subtle and complicated, as should become clear below, 
particularly in §3.1. 

The contraction of a neutron star due to spin-down also causes adiabatic heating. 
However, the adiabatic temperature change produced by a given density change A/9 is 
AT ~ (A/?/ /3)T, much smaller than the change from the process proposed above, for which 
AT ~ (A/?/ /3)/in, where /in is the chemical potential (or Fermi energy) of the neutrons (not 
including rest mass). Even smaller are the effects of the contraction due to cooling (Baym 
1981). 

This paper is organized as follows. In §2, I briefly discuss the equilibrium properties of 
a rotating star composed of neutrons, protons, and electrons, and their dependence on the 
rotation rate, O. In §3, I attempt a description of the chemical and thermal evolution of the 
star as O changes with time. In §3.1., I outline the (possibly very complicated) interplay 

^ The reaction rates used by Finzi & Wolf (1968) are not correct (Haensel 1992). Lai 
(1994) incorrectly distinguishes between "bulk viscosity" and "damping by neutrino emis- 
sion," which in fact correspond to the same physical process. 
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between reactions and diffusive processes in stars with nontrivial internal structure, such 
as superfluid and non-superfluid regions, or regions with substantiaUy different reaction 
rates. In the subsequent sections, the spatial structure of the star is ignored, except for 
the use of a standard formula relating the internal and surface temperatures. Section 3.2 
presents a very simplified model for the chemical and thermal evolution of a non-superfluid 
star, which reduces to two coupled flrst-order differential equations for the typical values of 
the departure from chemical equilibrium Sfx and the internal temperature T. The needed 
total neutrino and antineutrino emissivities and differential reaction rates are given in the 
Appendix. Section 3.3 contains the results of the numerical integration of the differential 
equations for stars spinning down due to magnetic dipole braking with different surface 
magnetic fleld strengths, and a qualitative analysis of these results. Changes to be expected 
with the introduction of pairing (superfluidity) are discussed in §3.4. In §4, I comment 
on the shortcomings of the present model, suggest possible improvements, compare the 
proposed mechanism to previously studied heating processes, and make a preliminary 
comparison of my results with recent observations. 



2. EQUILIBRIUM STATE OF A ROTATING STAR 

I consider a neutron star that is spinning down (or up) at a rate slow enough that 
it can always be considered to be in hydrostatic equilibrium and in uniform (rather than 
differential) rotation. Because of the high degeneracy of the neutron star matter, its 
cooling time and its relaxation time towards chemical equilibrium are much longer than 
its dynamical timescale, and also may be longer than the characteristic spin-down time. 
For simplicity, I use Newtonian dynamics throughout this paper. 

In hydrostatic equilibrium, 

= --Vp-V(^, (2) 
P 

where p is the mass density, p is the pressure, <f) = if) — (l/2)0^r^ sin^ 6 is the effective 
gravitational potential (including the centrifugal potential), ij) is the actual gravitational 
potential, O is the rotation rate, r is the radial coordinate, and 6 is the colatitude. Taking 
the curl of eq. (2), one sees that the surfaces of constant /?, constant p, and constant cf) 
coincide. It is convenient to label these surfaces by the amount of mass M enclosed by 
each of them, as S{M), and in this way be able to identify surfaces at different rotation 
rates, although the values of /?, p, and cf) on S{M) change as the star spins down. 

Since the matter in a neutron star core is composed of neutrons, protons, electrons, 
and other particles (muons, perhaps kaons, etc.), its equation of state must in general be 
written as p = p(/9, ) where the relative concentrations of different particle 

species. In this paper, I consider only neutrons, protons, and electrons. In this case, the 
proton number density and the electron number density are required to be equal in 
order to preserve charge neutrality. Thus, a single "composition parameter" x = rip/n is 
enough to characterize the chemical state of the matter (n = rij^ -\- is the total baryon 
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number density, and rij^ is the neutron number density). Since p and p are constant on 
each surface S{M), this must also be true for x. (If additional particles are present, more 
"composition parameters" become necessary, and these need not individually be constant 
on these surfaces in order for p and p to be constant simultaneously.) 

The chemical equilibrium state of the matter at a given baryon density, n, and given 
entropy per baryon, s, is found by minimizing the total energy per baryon, E{n,x,s) 
(including the contribution from the electrons), with respect to x, i.e., by demanding that 

dE 

— = = /ip + /ie - /in = 0, (3) 

where /in, /ip, and p^ are the chemical potentials of the neutrons, protons, and electrons, 
respectively. The equilibrium value of x (xeq) is a function of n. 

A change in the stellar rotation rate changes the centrifugal force, and therefore the 
baryon density n on a given surface iS(M), 

An aOAO , , 
~ • (4) 

n Gp, ^' 

for <C Gpc- Here, the symbol A denotes the Lagrangian variation (i.e., the variation 
on a given surface S{M)) of the quantity following it, G is the gravitational constant, 
Pc is the central density of the star, and a is a positive, dimensionless number of order 
unity that depends only slightly on radius within the core of a neutron star. In the center 
of a star with the analytically solvable polytropic equation of state p (x p^ (x (which 
falls into the range of theoretically allowed equations of state of neutron star matter), 
a = (vr^ — 3)/37r fii 0.73. Thus, the equilibrium value of x varies as 

dXeq nEnx OAO 

Axeq = — — An = a— — . (5) 

an Exx Gpc 

Here and in what follows, italicized subscripts on E are used to denote partial derivatives 
of the energy per baryon (e.g., E^x = E / dndx) to be evaluated in the cold, chemical 
equilibrium state (T = Sp = 0). 



3. CHEMICAL AND THERMAL EVOLUTION 

3.1. Reactions vs. diffusion 

Two effects make the actual value of x on S{M) change as O varies. First, since Xeq 
changes, the matter is out of chemical equilibrium, and this induces a difference V between 
the rates (per unit volume) of neutron beta decays and inverse beta decays that goes in 
the direction of restoring the equilibrium. Second, the total chemical potentials (including 
internal and external contributions, cf. Kittel & Kroemer 1980) of the individual particle 
species typically become different on different surfaces iS(M), thus deviating from diffusive 
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equilibrium. If the timescale for diffusion of particles is shorter than the timescale for 
relaxation to chemical equilibrium, significant diffusion occurs and the relaxation to the 
chemical equilibrium state is a global process involving the whole stellar core. In the 
opposite case, relaxation proceeds locally, without significant interactions among different 
parts of the star.^ Both regimes may occur in a real neutron star, depending on its history 
and on the uncertain properties of neutron star matter. The first is more likely if a 
superfluid is present (which inhibits reactions by reducing the available phase space, and 
probably also facilitates relative motions of the superfluid with respect to other particles), 
if the temperature is low, or if only modifled Urea processes are allowed. The opposite 
conditions tend to favor the second regime. 

For a given function r(n,(5/i,T), the flrst regime is more efficient in restoring equi- 
librium, since most reactions can occur in those regions where F is larger. For example, 
assume direct Urea reactions are allowed in the inner core region of the star but not in 
the outer core. If the diffusion time is short, Sfx is always similar in both regions, with 
particles diffusing from one to the other. Reactions occur almost exclusively in the inner 
region, and most of the heat is generated there. If the diffusion time is long, the inner 
core always stays closer to chemical equilibrium than the outer core. Reactions occur at 
comparable rates in both regions, and most of the heat generation and energy emission (in 
the form of neutrinos) occurs in the outer region. 

Gnedin and Yakovlev (1994) have proposed that a temperature gradient in the core 
of a neutron star may generate convective motions that would further complicate the 
thermal evolution. However, these authors do not take the strong stabilizing effect of the 
composition gradient into account, which wins by a large factor over the destabilizing effect 
of the temperature gradient (Pethick 1991; Reisenegger & Goldreich 1992), except possibly 
during the flrst few seconds after the supernova explosion (Thompson & Duncan 1993). 
The existence of a composition gradient also argues strongly against convection driven by 
a divergence of the magnetic stress tensor, which has been proposed by Urpin and Ray 
(1994). Radial motions of matter must be accompanied either by ambipolar diffusion (i.e., 
relative motions of different particle species) with relative velocities comparable to the bulk 
displacements, or by fast reactions that keep the moving fluid always close to its current 
chemical equilibrium state (Pethick 1991; Goldreich & Reisenegger 1992). 

The model studied in the following sections ignores all these complications. It is 
unclear at this point how their inclusion would affect the predicted surface temperature 
of the star. This uncertainty should be kept in mind when analyzing the results of the 
calculations reported below. If possible, it should be addressed by future, more reflned 
calculations. 

3.2. A simple model 

In this section, I present a simple model for the thermal evolution of a neutron star in 
the presence of internal chemical heating due to spin-down. I take the matter to consist of 

^ These two regimes are also discussed more quantitatively, in the context of magnetic 
fleld evolution in neutron stars, by Goldreich & Reisenegger (1992). 
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non-superfluid neutrons, protons, and electrons, and assume that the chemical and thermal 
state of the star can be described by just two variables, which are assumed to be indepen- 
dent of position within the star: the departure from chemical equilibrium defined in 
eq. (3), and the internal temperature T. When needed, I calculate the effective surface 
temperature Tg from the internal temperature by the relation derived by Gudmundsson, 
Pethick, and Epstein (1982). The model includes the variation of the chemical equilibrium 
state due to spin-down, the relaxation to chemical equilibrium, heating or cooling due to 
neutrino-emitting reactions, and cooling due to emission of photons from the surface. 

For cool, degenerate matter close to chemical equilibrium (in the sense that the total 
energy differs little from its value for perfectly cold matter in chemical equilibrium at the 
same baryon density), the total energy per baryon can be written as 

E{n, X, s) = Eo{n) + Echin, x) + Ethin, s), (6) 
where Eq is the energy in the cold, chemical equilibrium state of baryon density n, 

Ech = IeU.x - xe^f = lE^^'i^f^f (7) 
is the "chemical" energy (due to the departure from chemical equilibrium), and 

^th = Ie,,s' = ^E-'T' (8) 

is the thermal energy. The second equality in eqs. (7) and (8) follows from the ther- 
modynamic identities (dE / ds)n,x = T and (dE / dx)n,s = Sfx. The time evolution of the 
chemical energy is given by 

Ech = Exx{x - Xeq){x - igq) = E~l 8 jlS jl 

y Exx Gpc n j ' 

In the last expression, the first term in parentheses accounts for the change in the chemical 
equilibrium state due to spin-down (or spin-up), and the second term for the change in the 
actual chemical state due to reactions. Similarly, the thermal evolution is determined by 



E^. = E-}TT 



= r.^^ 11 = E^, (10) 

n n 

where the three terms in the last expression are, in the same order, the chemical energy 
released by the reactions, the energy radiated in neutrinos and antineutrinos, and the 
energy radiated in photons from the surface. Expressions for VSiJ. and e are given in the 
Appendix. For the last term, I make the rough estimate 

E^ ^ -^L^ ^ 4.5 X IQ-^^T^^ eTgs-\ (11) 
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where mn and are the masses of a neutron and of the star, and the surface photon 
luminosity, L-y, is given in terms of the interior temperature, T = 10^ Tg K, by the relation 
of Gudmundsson et al. (1982), claimed to hold to better than 1.5% in the range 0.16 < 
Tg < 10. For illustration purposes, I assume that it also holds for temperatures beyond 
the lower limit of this range. Since the lower limit is already lower than the typical current 
detection threshold, this should not cause serious problems. 

For Eas^ I use the value corresponding to a mixture of non-relativistic neutrons and 
protons and extremely relativistic electrons, all non-interacting and completely degener- 
ate: 



e: 



2/3 /k^ 



(1 



-(Svr^n) 



l/3™2/3 
■^eq 



(12) 



In practice, interactions among nucleons are expected to reduce their effective mass relative 
to mn, but this effect is uncertain and unlikely to change the result by a large factor, so it 
is ignored here. 

In order to estimate Enx and -E^x, I use the customary (and very good) approximation 
(see Prakash, Ainsworth, & Lattimer 1988 for references) 

E^h(n, x) = E{n) + S{n){l - 2xf + ^nc{?,^'^nfl^x^l^, (13) 

where E{n) = —S{n){l — 2xeq)'^ — (3 / 4)?ic{37r'^ riY xtq^ , S{n) is the so-called "nuclear 
symmetry energy," and the last term is the contribution of the electrons. From eqs. (6) 
and (13), I obtain 

Enx = -(q-l) nc{37r'x,q)'/'n-'/\ (14) 



au' 



d 



/o 2 \l/3 -2/3-'-+4^eq /^P-x 

Exx = -nc(3T: n) ' x^„ ' , 15 

with q = din S I dlnn, and Xeq given by the condition E^ix = x^q) = 0, or equivalently 



_ 45(n) 

l-2xeq ;ic(37r2n)i/3 • ^ > 

In what follows, I take a = 0.73, pc = 3m^nf) = 8.0 X lO'^'^gcm"^, n = 2no, and 
Xeq = 1/9, the latter corresponding to the threshold for direct Urea processes to be allowed 
(Lattimer et al. 1991). The value of q (see eq. [14]) is particularly uncertain. It can not 
even be decided with certainty whether q is larger or smaller than 1/3. This determines 
the sign of Enx and dxeq/dn^ and therefore of the departure from chemical equilibrium. I 
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assume q = 2/3, the value for non-interacting, degenerate nucleons. In addition, I follow 
Taylor, Manchester, and Lyne (1993) in defining the spin-down power as 



V = 10'*^P36 erg s-i = -/OO, (17) 

with an assumed moment of inertia / = IC^^gcm^. Table 4 of Taylor et al. (1993) 
shows the observed values of V for a large number of radio pulsars. I also write Sfx = 
irk X 10^(^8 K ^ 27(^8 keV, and u = Sg/Tg = 8 ji / {T^kT). With these conventions, and the 
rates given in the Appendix, the evolution equations (9) and (10) become 



^8 = (-1.0^36 - 3.8 X IQ'SgTf—^^^^^^^ Myr-^ 

• / „ 7_,21m6 + lOSu^ - 357^2 - 457 _i 2 \ . . 1 

T8 = f 1.7 X WtI — QmTl^ ] Myr-\ 

if direct Urea reactions are active, and 



(18) 



A f in-p , ^6 + 105u^ + 945u^ + 1835 ^ 
Ss = { -1.0^36 - 2.5SsT^ ^^^^ 1 Myr , 

• , _7 15m^ +420m*5 + 1890m^ - 7340^2 - 11513 _i 2 \ . . 1 

Ts= (0.14T; 0.92T^ ') Myr-\ 



(19) 



if only modified Urea reactions are allowed. The ratios of the coefficients of different 
powers of u in each equation are exact in the limit of vanishing Sfx and T (see Appendix 
and Haensel 1992), but all other coefficients in these equations are only rough estimates. 
Note that a modest departure from chemical equilibrium, \u\ < 1.5, increases the cooling 
rate, and only larger departures from chemical equilibrium reduce the cooling rate or may 
even lead to net heating. For both direct and modified Urea processes, the cooling rate is 
enhanced by a maximum factor ^1.5 for \u\ ^ 1.1. 



3.3. Results 

Neutron stars are born hot (T ~ 10"^"^ K), but with rotation rates and magnetic fields 
not enormously larger than the observed ones. Thus, in the earliest phase of the evolu- 
tion, spin-down heating and photon cooling are negligible, and only Urea reactions are 
important. Under these assumptions, 

,'3-1^^10 (20) 

for both direct and modified Urea reactions. This result implies that \Sfj,\ decreases much 
faster than T, and soon only the lowest-order powers of u contribute to the evolution. In 
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this case, T oc t~^/'^ for direct Urea reactions, T oc t~^/^ for modified Urea reactions, and 
oc t~l^ for both. Thus, as has been argued before (e.g., MaxweU 1979), but (to my 
knowledge) never shown exphcitly, a primordial departure from chemical equilibrium in a 
nonrotating star decays on a timescale substantially shorter than the cooling time, unless 
Urea reactions are not the dominant cooling mechanism, and the latter does not alter the 
chemical composition. 

In a star that is spinning down, increases again after the temperature (and hence 
the reaction rate) is low enough. If the spin-down power can be approximated as a constant 
(which is the case at early times for a pulsar with magnetic dipole braking), one has in 
this stage 



^ (21) 

1 + /3' Myr ' ^ ^ 



with /3' = /3/4 fii 6 for direct Urea reactions, and /3' = /3/6 ~ 3 for modified Urea reactions, 
according to equations (18) and (19). increases in this way until it is of order A;T, at 
which point higher order powers of u become important in the expressions for the reaction 
rates, and the conversion of chemical energy into thermal energy is speeded up. Only 
after this time does the finite departure from chemical equilibrium substantially affect the 
thermal evolution. 

If V varies slowly over the timescales required to cool the stellar core and achieve 
chemical equilibrium, then a quasi-equilibrium state is reached in which the terms on 
the left-hand sides of equations (9) and (10) are negligible compared to individual terms 
on the right-hand side, which nearly cancel each other. In this state, the departure from 
chemical equilibrium and the core and surface temperatures are independent of the previous 
evolution of the star, and are completely determined by the current value of "P, which is 
directly measurable from pulsar timing, and by the neutron star model (mass, composition, 
and equation of state). 

Fig. 1 shows numerical integrations of equations (18) and (19) together with the corre- 
sponding quasi-equilibrium solutions. I assume magnetic dipole braking, PP ^ 10~^^i?^, 
with P in seconds and the surface magnetic field strength B in gauss, and I assume a 
fast initial rotation. Pi = 1 ms, which maximizes the heating due to spin-down at early 
times, t < 10^^(Pi [ s]/i?[ G])^ s. The three panels correspond to three different (constant) 
magnetic field strengths, (a) B = 10^^ G, (b) B = 10^° G, and (c) B = 10^ G. For the 
strongest field (10"^^ G, representative of a typical, "classical" radio pulsar), the spin-down 
occurs so rapidly that the heating is only substantial at early times, when the temperature 
is still quite high. Thus, the only observable effect is a small bump in the cooling curve, 
which would be even smaller for a (perhaps more realistic) lower initial rotation period. 
The heating effect becomes more interesting for the intermediate field (10"^" G). Here, the 
evolutionary curves follow the quasi-equilibrium curves very closely after Sfx has grown 
enough for the quasi-equilibrium state to be reached for the first time. In the photon- 
cooling regime, some deviation is seen, because the reactions are too slow to follow the 
rapidly dropping equilibrium solutions. For the lowest field (10^ G, representative of a mil- 
lisecond pulsar), the heating is negligible until around the time when the photon-cooling 



10 



regime is reached. After this point, however, both the departure from chemical equihb- 
rium and the temperature foUow the quasi-equihbrium curves, which remain quite high 
for about a Hubble time. Of course, the "standard model" for the history of millisecond 
pulsars (see Bhattacharya & van den Heuvel 1991 for a review) is not as uneventful as the 
evolution of my 10^ G neutron star. The magnetic field may have decayed or the star may 
have been heated by accretion. Still, unless these events were very recent, the star should 
have evolved to the quasi-equilibrium state, and its luminosity should be predictable from 
its spin-down power. 

Fig. 2 shows the evolution of the effective surface temperature for B = 10^^ G (dotted), 
10^^ G (short-dashed), 10^° G (long-dashed), 10^ G (short-dashed-dotted), and 10^ G (long- 
dashed-dotted), and in the absence of spin-down heating (solid), for direct (lower curves) 
and modified (upper curves) Urea reactions. Again, one sees that the effect of spin-down 
heating can be substantial at late times for stars with weak fields. The effect is quite weak 
for stars with strong fields, and is further reduced if the initial period is longer than the 
assumed 1 ms. 

3.4- Speculations on superfluid effects 

So far, there is no published calculation of the reaction rates in superfluid (i.e., BCS- 
paired) neutron star matter out of chemical equilibrium. Thus, it is not yet possible to 
compute the coupled thermal and chemical evolution of a superfluid neutron star, not 
even in a very simplifled way such as that presented in the previous sections for a normal 
neutron star. In the following paragraph, I make educated guesses about the qualitative 
changes expected in the evolution if pairing effects were taken into account. 

The thermal energy at a given temperature (or, equivalently, the heat capacity) of 
superfluid particles is substantially reduced if the temperature is much smaller than the 
transition temperature Tc- The chemical energy, however, should be nearly unaffected 
by the pairing. The reactions involving superfluid particles are substantially slowed if 
T <C Tc and <C kTc- The reduction in the reaction rates makes \Sfj,\ rise faster, and 
allows a larger fraction of the rotational energy to be stored as chemical energy. For the 
same reason, photon cooling becomes important at a higher temperature, and the ratio 
\u\ = \Sfj,\/TrkT is larger. The larger chemical energy is stored for a longer time, allowing 
heating to take place until later in the evolution of the star. In the quasi-equilibrium 
regime, the reduced reaction rate increases the photon luminosity for a given spin-down 
power, just as in normal matter the equilibrium photon luminosity for the modifled Urea 
processes is substantially higher than that for the much faster direct Urea processes. 

These arguments suggest that the effects of spin-down heating are more dramatic in 
a superfluid neutron star than in the normal star studied more carefully in this paper. 

4. DISCUSSION 

The quantitative results presented in the previous section are based on a very sim- 
plifled (perhaps oversimplifled) model that ignores many features to be expected in a real 
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neutron star. These include a non-trivial spatial structure, where different layers can have 
different properties, and where diffusion from one layer to another (or its absence) can 
have an important effect. Related to this is the finite timescale for thermal conduction, 
which is known to modify the early stages of neutron star cooling (Nomoto & Tsuruta 
1987). For many of the physical parameters affecting the evolution, I took representative 
values which may well be off by a factor ~ 2 or more. I also considered only two cooling 
mechanisms (neutrino emission by direct or modified Urea reactions, and surface photon 
emission), neglecting several others (see Nomoto & Tsuruta 1987). Finally, superfluidity 
is not included in the model. Its effects are discussed briefly in §3.4. 

Two main pieces of work would be needed in order to produce models that could 
satisfactorily explore the thermal evolution of different model neutron stars. First, in 
order to include neutron superfluidity and proton superconductivity, it will be necessary to 
extend the reaction rate calculations to superfluid neutron star matter at flnite temperature 
and away from chemical equilibrium, possibly building on recent work by Haensel (1992) on 
normal matter at flnite temperature and away from chemical equilibrium, and by Levenflsh 
& Yakovlev (1994b) on superfluid matter at flnite temperature and in chemical equilibrium. 
Second, one should use a code with spatial structure in at least one (radial) dimension, in 
order to take into account the diffusion of particles and heat and to make a quantitative 
study of some of the effects pointed out in §3.1. Of course, a rotating star is not spherically 
symmetric but, as I showed in §2, the variables determining the stellar structure are all 
constant on one set of surfaces, S{M). This probably makes a one-dimensional description 
practical. The code could include all known cooling mechanisms, and the effect of General 
Relativity, in order to produce predictions which could be compared directly and reliably 
to observations (to the extent that the latter are reliable as well; see below). 

In spite of its obvious shortcomings, the simple model shows the main qualitative 
effects of the conversion of rotational energy into chemical energy, and flnally into heat. 
Since the rate at which energy is converted at each of the two steps depends on the values 
of Sfi and T, it cannot simply be treated as a given source term in the thermal evolution. 
In this respect, it is very different from and somewhat more subtle than the other heating 
mechanisms studied so far and mentioned in the Introduction. However, it may be regarded 
as less speculative than the previously proposed heating mechanisms based on vortex creep 
or crust cracking, since it involves only physical processes and parameters which affect the 
thermal evolution even in the absence of heating, and not additional (and very uncertain) 
parameters such as the pinning strengths of vortex lines to the solid lattice of the crust and 
the yield strength of the latter. It is often attempted to determine the values of the latter 
by comparing theoretical models for the thermal evolution to observations (e.g., Cheng et 
al. 1992; Van Riper et al. 1994). Since all heating mechanisms are most noticeable at late 
times in the thermal evolution of a neutron star, these studies are subject to confusing 
different heating mechanisms, and their observational disentanglement does not appear 
feasible in the near future. 

It is clear from the numerical integrations (§3.3.) that the proposed new heating effect 
is only of small to moderate importance in the early (neutrino-dominated) era of the ther- 
mal evolution, but can become quite dramatic in later stages, particularly if the interior of 
the star is superfluid (see §3.4.). It should be interesting to compare the quasi-equilibrium 
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luminosities predicted for millisecond pulsars to observational data. Unfortunately, it is 
hard to extract the surface luminosity due to the heat flowing out of the stellar interior 
("cooling luminosity") from the observations. First, interstellar absorption modifles the 
spectra and observed luminosities and is hard to correct for, since the intrinsic spectra 
are not known (e.g., Becker & Triimper 1993; Banner, Kulkarni, & Thorsett 1994). Sec- 
ondly, the luminosity is usually (perhaps always) dominated by the additional emission 
coming from the star's magnetosphere or from external (magnetospheric) heating of the 
polar caps (see, e.g., Ogelman 1994 for a review). Based on "classical" pulsars, Seward 
& Wang (1988) and Ogelman (1994) have proposed empirical formulae for the magneto- 
spheric x-ray luminosity in terms of the spin-down power, which predict that the former 
should indeed dominate over the "cooling luminosity" by a large factor (Fig. 3). However, 
Banner et al. (1994) have pointed out that the luminosities of some millisecond pulsars 
(detections or upper limits) are much lower than predicted by the empirical relations of 
Seward & Wang (1988) and Ogelman (1994), and in fact differ signiflcantly among pulsars 
with very similar spin-down power. This gives hope that some millisecond pulsars may 
provide a "window" through which the "cooling luminosity" can be seen. 

Fig. 3 shows that, if direct Urea processes are permitted, the luminosity originating 
from the interior should be completely unobservable. This is not the case, however, if 
only modifled Urea processes are allowed, for which the predicted cooling luminosities 
come quite close to some of the ROSAT PSPC upper limits of Banner et al. (1994). Of 
course, in view of the deflciencies of the present model and of the uncertainties in the 
observations, this result can only be taken to indicate that the proposed effect may well 
be interesting, and that further theoretical and observational investigations are warranted. 
The former have already been discussed above, and among the latter the most interesting 
may be x-ray (or extreme ultraviolet) observations of the millisecond pulsar B1937-|-21, 
whose spin-down power is only a factor of two lower than that of B1821-24 (see Fig. 3), 
and which is probably closer (~ 3.6 kpc compared to ~ 5.1 kpc; Taylor et al. 1993, Banner 
et al. 1994). 

5. CONCLUSIONS 

The spin-down of a neutron star perturbs its internal chemical equilibrium state, 
leading to reactions and to internal heating. Its effect on the thermal history of stars 
with strong magnetic flelds (> 10"^"^ G) is not very large, but it becomes substantial for 
weak-fleld stars, especially at late times. Old neutron stars reach a quasi-equilibrium state 
in which heating and cooling balance each other, and (for a given neutron star model) the 
surface photon luminosity can be predicted from the spin-down power. If only modifled 
Urea reactions can occur, this luminosity lies in the observable range for at least some 
known millisecond pulsars. More detailed calculations, including the determination of 
reaction rates for superfluid matter out of chemical equilibrium, are needed in order to 
make precise predictions for speciflc stellar models. 
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the initial idea for this paper. I also thank Rudi Banner and Shri Kulkarni for useful 
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APPENDIX 

REACTION RATE AND NEUTRINO EMISSIVITY 



In the limit where kT and Sfx are much smaller than the Fermi energies of neutrons, 
protons, and electrons, and assuming that all particles are "normal" (non-superfluid), 
Haensel (1992) showed that the emissivities (energy per unit volume per unit time emitted 
in neutrinos/antineutrinos) can be written as^ 

e^iT,Sf,) = e,{T,-6^) = 0)^^^^|^, (22) 

and the emission rates (number of neutrinos/antineutrinos emitted per unit volume per 
unit time), as 

T^{T,Sn) = Ti,{T,-Sn) = — • (23) 

Here, e(T, Sfx) = ejj + eu denotes the total emissivity, including neutrinos and antineutrinos, 

/"^ , y^P(y-x) , , 

F{x) ^ / dy-^i-^ L_ 24 

Jo l + exp(y-a;) 

and 

For the direct Urea processes, P{y) = -Pd(y) = ^^"^ + y^7 ^-nd for the modified Urea processes, 
P(y) = P^(y) = 97r4 + 107r2y2+y4_ 

For the present work, only the total emissivity e and the net reaction rate F = Fj^ — Fj^ 
are needed. These are proportional to 



F+{X) ^ F{X) + F{-X) = r dyi^ll+^piinyl _ r ^y^y _ ^^p^y^^ ^^Q^ 

Jo 1 + Jo 



^ The present sign convention for Sji is that used by Reisenegger & Goldreich (1992) 
and Goldreich & Reisenegger (1992), and opposite to that of Sawyer (1989) and Haensel 
(1992). 
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and 



G-{x) = G{x) - G{-x) = 4x / t/yf^ + / dy{y - xfP{y), (27) 

where I have used the fact that P{y) = P{—y). Using the exphcit forms for P{y) given 
above, and the identity 



= ,08, 

l + ey V 22'^- V An ^ ' 

(Gradshteyn & Ryzhik 1980, eq. 3.411.4), where Bn are the BernouUi numbers (i?2 = 1/6, 

i?4 = —1/30, i?6 = 1/42, i?8 = —1/30), one can write F^{x) and G~{x) exphcitly as 
polynomials. Defining u = 8 [j. I {irkT), I obtain 

^rr. c s ^rr. J 1 7 1 M ^ + 3 1 + 2 1 M ^ \ 

ed(T, 8^) = eaiT, 0) fl + — j (29) 



and 



Td{T, SfijSfi = ed(T, 0) — (30) 

457 



for the direct Urea processes, and 



em(T, Sn) = em(T, 0) ( 1 + ^j^j^ ) (31) 



and 



^ , , , , 14680m2 + 7560m^ + 8A0u^ + 24m« , , 

r^(T, ,5/i),5/i = e^(T, 0) (32) 

for the modified Urea processes. 

Throughout this paper, I use the (somewhat uncertain) equilibrium emissivities given 
by Haensel (1992; references to the original derivations and to critiques are given there). 



1/3 

ed(T, 0) ^ 4.3 X 10^1 ( ] erg cm'^ s'^ (33) 



no 

and 



1/3 

no 



e^(T, 0) ^ 3.5 X 10" ( ^ ) T« erg cm'^ s'^ , (34) 



where no = 0.16 fm ^ is the normal nuclear baryon density. 



15 



REFERENCES 



Alpar, M. A., Anderson, P. W., Pines, D., k Shaham, J. 1984, ApJ, 276, 325 

Bahcall, J. N. & Wolf, R. A. 1965, Phys. Rev., 140, B1452 

Baym, G. 1981, ApJ, 248, 767 

Becker, W., & Triimper, J. 1993, Nature, 365, 528 

Bhattacharya, D. & van den Heuvel, E. P. J. 1991, Phys. Rep., 203, 1 
Boguta, J. 1981, Phys. Lett. B, 106, 255 

Brown, G. E., Kubodera, K., Page, D., & Pizzochero, P. 1988, Phys. Rev. D, 37, 2042 

Cheng, K. S., Chau, W. Y., Zhang, J. L., & Chau, H. F. 1992, ApJ, 396, 135 

Chiu, H. Y. & Salpeter, E. E. 1964, Phys. Rev. Lett., 12, 413 

Cutler, C, Lindblom, L., & Splinter, R. J. 1990, ApJ, 363, 603 

Banner, R., Kulkarni, S. R., & Thorsett, S. E. 1994, ApJ, in press 

Finzi, A. & Wolf, R. A. 1968, ApJ, 153, 835 

Gamov, G. & Schoenberg, M. 1941, Phys. Rev., 59, 539 

Gnedin, O. Y. & Yakovlev, D. G. 1994, Astron. Lett., 20, 40 

Goldreich, P. & Reisenegger, A. 1992, ApJ, 395, 250 

Gourgoulhon, E. & Haensel, P. 1993, A&A, 271, 187 

Gradshteyn, L S., & Ryzhik, L M. 1980, Table of Integrals, Series, and Products, fourth 
edition (New York: Academic Press) 

Gudmundsson, E. H., Pethick, C. J., & Epstein, R. L 1982, ApJ, 259, L19 

Haensel, P. 1992, A&A, 262, 131 

Haensel, P. & Schaeffer, R. 1992, Phys. Rev. D, 45, 4708 

Kittel, C. & Kroemer H. 1980, Thermal Physics, second edition (New York: Freeman), 
chapter 5 

Lai, D. 1994, MNRAS, in press 

Landau, L. D., & Lifshitz, E. M. 1959, Fluid Mechanics (Oxford: Pergamon), §§15 & 78 

Lattimer, J. M., Pethick, C. J., Prakash, M., Haensel, P., 1991, Phys. Rev. Lett. 66, 2701 

Levenfish, K. P. & Yakovlev, D. G. 1994a, Astron. Rep., 38, 247 

Levenfish, K. P. & Yakovlev, D. G. 1994b, Astron. Lett., 20, 43 

MaxweU, O. V. 1979, ApJ, 231, 201 

Nomoto, K. & Tsuruta, S. 1987, ApJ, 312, 711 

Ogelman, H. 1994, in Lives of the Neutron Stars, NATO ASI series, in press 

Pethick, C. J. 1991, in The Structure and Evolution of Neutron Stars, ed. D. Pines, R. 
Tamagaki, & S. Tsuruta (Redwood City: Addison- Wesley), 115 



16 



Pethick, C. J. 1992, Rev. Mod. Phys., 64, 1133 

Prakash, M., Ainsworth, T. L., & Lattimer, J. M. 1988, Phys. Rev. Lett., 61, 2518 
Reisenegger, A. & Goldreich, P. 1992, ApJ, 395, 240 
Sakurai, T. 1969, J. Fluid Mech., 37, 689 

Sawyer, R. F. 1989, Phys. Rev. D, 39, 3804; erratum 1989, Phys. Rev. D, 40, 4201 
Seward, F. D. & Wang, Z.-R. 1988, ApJ, 332, 199 
Shibazaki, N. & Lamb, F. K. 1989, ApJ, 346, 808 

Taylor, J. H., Manchester, R. N., & Lyne, A. G. 1993, ApJSS, 88, 529 
Thompson, C. & Duncan, R. C. 1993, ApJ, 408, 194 

Umeda, H., Shibazaki, N., Nomoto, K., & Tsuruta, S. 1993, ApJ, 408, 186 
Urpin, V. A. & Ray, A. 1994, MNRAS, 267, 1000 

Van Riper, K. A., Link, B., k Epstein, R. L 1994, preprint (astro-ph 9404060) 
Walin, G. 1969, J. Fluid Mech., 36, 289 
Wolf, R. A. 1966, ApJ, 145, 834 



17 



FIGURE LEGENDS: 

Fig. 1. — Thermal and chemical evolution of neutron stars. The variables Sfx/irk (a mea- 
sure for the departure from chemical equilibrium; short-dashed lines) and T (the internal 
temperature of the star; long-dashed lines) are plotted logarithmically as functions of time, 
with direct (lower curves) and modified (upper curves) Urea reactions. The dotted curves 
are the quasi-equilibrium solutions obtained at each point by requiring that the right-hand 
sides of the evolutionary equations (18) and (19) vanish. For reference, the temperature in 
the absence of heating is also shown (solid line). The three panels correspond to magnetic 
dipole braking with different magnetic field strengths: (a) B = 10^^ G; (b) B = 10^^ G; 
(c) B = 10^ G. The initial spin period is taken to be 1 ms. 

Fig. 2. — Effective surface temperature as a function of time for stars with direct (lower 
curves) and modified (upper curves) Urea reactions, with no heating (solid lines) or spin- 
down heating with magnetic field strengths B = 10^^ G (dotted), 10"^"^ G (short-dashed), 
10"^" G (long-dashed), 10^ G (short-dashed-dotted), and 10^ G (long-dashed-dotted). The 
initial spin period is taken to be 1 ms. 

Fig. 3. — Surface photon luminosity of a neutron star as a function of its spin-down 
power. The solid line shows the bolometric luminosity as predicted by the model described 
in this paper for a star whose core consists of normal neutrons, protons, and electrons, 
and in which only modified Urea reactions are allowed. The long-dashed line shows the 
corresponding predicted luminosity in the ROSAT PSPC range (0.1 to 2.4 keV), assuming 
a blackbody spectrum. The short-dashed curve is a rough prediction for the bolometric 
luminosity with direct Urea reactions. For comparison, I show the empirical curves of 
Seward & Wang (1988; dot-dashed) and Ogelman (1994; dotted) for the magnetospheric 
emission, and the ROSAT PSPC observations of millisecond pulsars from Banner et al. 
(1994). Triangles indicate upper limits, asterisks detections. The filled triangle is the 
upper limit and likely detection of PSR B1821-24, corrected for absorption. 
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